#!/usr/bin/perl
# 2013/11/10 Categorize significant peaks in several features specified by BEDs

use strict;
use warnings;
use IO::File;
use Getopt::Std;
use Carp qw(croak);
$| = 0;


use File::Basename;

if($ARGV[0] eq '--help'){
	die "Usage : $0 [siginificant peaks in BED] [several target BED files]\n";
}

# どれくらい外に広げるか(片側の広げる距離）
my $SEARCH_OUTSIDE = 150;

my $FILE_significant = shift @ARGV;


#---------------------------------------
# BED fileを読み込む
#---------------------------------------
my %gene;
my %BED_NAMES;
my $NUMBER_FILES = 0;
foreach my $file(@ARGV){
	my ($name, $dir, $ext) = &fileparse($file, '\..*');

	$NUMBER_FILES++;
	$BED_NAMES{$name} = $NUMBER_FILES;

	my $fh_in = IO::File->new($file) or die "cannot open $file: $!";
	while($_ = $fh_in->getline()){
		s/\r?\n//;
		my ($chr, $start, $end, @others) = split /\t/;
		for(my $i = $start; $i <= $end; $i++){
			push @{$gene{$chr}{$i}}, $name;
		}
	}
	$fh_in->close();
}


my %category;
my $fh_sig = IO::File->new($FILE_significant) or die "cannot open $FILE_significant: $!";
while($_ = $fh_sig->getline()){
	s/\r?\n//;
	my ($chr, $start, $end, @others) = split /\t/;

	### Anotationと関連付ける
	my $anotation = &connectAnnotation($chr, $start - $SEARCH_OUTSIDE, $end + $SEARCH_OUTSIDE);

	$category{$anotation}++;
}
$fh_sig->close();

print "category\tnumber\n";
foreach my $key(sort {$category{$b} <=> $category{$a}} keys %category){
	print "$key\t$category{$key}\n";
}



sub connectAnnotation{
	my ($chr, $start, $end) = @_;
	my %candidate = ('NA' => 0);
	for(my $i = $start; $i <= $end; $i++){
		if(exists $gene{$chr}{$i}){
			foreach my $t(@{$gene{$chr}{$i}}){
				$candidate{$t}++;
			}
		}
	}
	# 最もオーバーラップしていた遺伝子を取得する
	foreach my $key(sort {$candidate{$b} <=> $candidate{$a}} keys %candidate){
		return $key;
	}
}